Influence of the bins

Author

Lisa Nicvert

Published

June 29, 2023

This code infers the parameters of a Hawkes process on the same data as in folder 05_example_real_data but with different bins widths.

Load libraries

# Load functions from main folder
library(camtrapHawkes)

packages <- c("here", "lubridate", "UnitEvents",
              "ggplot2", "patchwork", "dplyr")
base::lapply(packages, require)
here::i_am("analyses/07_bins_width/bins_width.qmd")

Parameters

# --- Inference parameters
# Bins width
dhours_vec <- c(3, 9)

ndays <- 1.5
gamma <- 0.5
scale <- 10e7
estimator <- "BVL"

# --- Silhouettes path
silhouettes_path <- here("data/species_silhouettes/")
silhouettes_labs <- c(impala = paste0("<img src='", silhouettes_path, 
                          "impala.png'", " height='20'/>", 
                          "<br>*impala*"),
                      kudu = paste0("<img src='", silhouettes_path, 
                                    "kudu.png'", " height='20'/>", 
                                    "<br>*kudu*"),
                      lion = paste0("<img src='", silhouettes_path, 
                                    "lion.png'", " height='20'/>", 
                                    "<br>*lion*"),
                      wildebeestblue = paste0("<img src='", silhouettes_path, 
                                              "wildebeest.png'", " height='20'/>",
                                              "<br>*wildebeest*"),
                      zebraburchells = paste0("<img src='", silhouettes_path, 
                                              "zebra.png'", " height='20'/>",
                                              "<br>*zebra*"))

# --- Figures
figures_path <- here("figures/07_bins_width")

Read data

dat <- read.csv(here("outputs/05_example_real_data/cleaned_data.csv"))

# Format datetimes
fmt <- "%Y-%m-%d %H:%M:%S"
dat$datetime <- as.POSIXct(dat$datetime, format = fmt)
dat$datetime <- force_tz(dat$datetime, tz = "Etc/GMT-2")

Infer several models

# --- get alphabetical species order 
spp_names <- sort(unique(dat$snapshotName))

modlist <- vector("list", 
                  length(dhours_vec))
names(modlist) <- dhours_vec

for (i in 1:length(dhours_vec)) {
  # --- Get k and delta
  delta <- dhours_vec[i]/24 # delta in days
  k <- ceiling(ndays/delta) # k from delta and dmax
  
  # --- Infer
  inf <- infer(dat, 
               k = k, delta = delta,
               gamma = gamma, scale = scale)
  # If the error:
  # Erreur dans NeuroCorr_SetParallel(threads) : 
  #   impossible de trouver la fonction "NeuroCorr_SetParallel"
  # appears, restart R, re-load packages and re-run
  
  # --- Format inferred model
  u_inf <- unpack_inf(inf, 
                      species_names = spp_names, 
                      k = k, delta = delta)
  
  
  df_ue <- ue_model_to_df(u_inf[[estimator]])
  
  modlist[[i]] <- df_ue
}

Plot inferred interactions and rates

Now let’s plot the graphs:

# --- Initialize glist
glist_inter <- vector("list", length(modlist))
names(glist_inter) <- names(modlist)

for (i in 1:length(modlist)) {
  # --- Plot
  mod <- modlist[[i]]
  bin <- names(modlist)[i]
  
  g <- plot_interactions(mod,
                         scale = "hours", 
                         textsize = 15,
                         timestep = floor((max(mod$time)*24)/3),
                         h = 0,
                         separate_self = TRUE,
                         silhouettes = silhouettes_labs)
  
  glist_inter[[i]] <- wrap_elements(g)
}
# --- Initialize glist_bg
glist_bg <- vector("list", length(modlist))
names(glist_bg) <- names(modlist)

for (i in 1:length(modlist)) {
  # --- Plot
  mod <- modlist[[i]]
  bin <- names(modlist)[i]
  
  g <- plot_background_rate(ue_df = mod,
                     textsize = 15,
                     write_label = TRUE,
                     silhouettes = silhouettes_labs,
                     title = "Background rates") +
    ylim(0, 0.25)
  
  glist_bg[[i]] <- wrap_elements(g)

}

Comparison of different bins widths

Here are the parameters inferred with 6 hours bins in 05_example_real_data:

Original interactions

Original background rates

Below are the re-infered values:

plot_names <- as.character(dhours_vec)
for (pn in plot_names) {
  wrap_plots(c(glist_inter[pn],
             glist_bg[pn]), 
           nrow = 2, byrow = TRUE,
           heights = c(8, 4))
  
  # Export plot
  ggsave(file.path(figures_path, paste0("binwidth_", pn, ".jpeg")),
         bg = "white",
         width = 9, height = 14,
         dpi = 300)
}